Journal  of  Power  Sources  193  (2009)  258-268 


ELSEVIER 


Contents  lists  available  at  ScienceDirect 

Journal  of  Power  Sources 

journal  homepage:  www.elsevier.com/locate/jpowsour 


pri 

SllbiidtS 


Operational  optimization  and  real-time  control  of  fuel-cell  systems 

J.  Hasikos,  H.  Sarimveis* *,  P.L.  Zervas,  N.C.  Markatos 

National  Technical  University  of  Athens,  School  of  Chemical  Engineering,  9  Heroon  Polytechniou  Str.,  GR-157  80  Athens,  Greece 


ARTICLE  INFO 


ABSTRACT 


Article  history: 

Received  10  October  2008 

Received  in  revised  form  7  January  2009 

Accepted  8  January  2009 

Available  online  29  January  2009 


Keywords: 

Fuel  cells 
Hydrogen 

Model  predictive  control 
Optimization 
Meta-modeling 
Neural  networks 


Fuel  cells  is  a  rapidly  evolving  technology  with  applications  in  many  industries  including  transportation, 
and  both  portable  and  stationary  power  generation.  The  viability,  efficiency  and  robustness  of  fuel-cell 
systems  depend  strongly  on  optimization  and  control  of  their  operation.  This  paper  presents  the  devel¬ 
opment  of  an  integrated  optimization  and  control  tool  for  Proton  Exchange  Membrane  Fuel-Cell  (PEMFC) 
systems.  Using  a  detailed  simulation  model,  a  database  is  generated  first,  which  contains  steady-state 
values  of  the  manipulated  and  controlled  variables  over  the  full  operational  range  of  the  fuel-cell  sys¬ 
tem.  In  a  second  step,  the  database  is  utilized  for  producing  Radial  Basis  Function  (RBF)  neural  network 
“meta-models”.  In  the  third  step,  a  Non-Linear  Programming  Problem  (NLP)  is  formulated,  that  takes 
into  account  the  constraints  and  limitations  of  the  system  and  minimizes  the  consumption  of  hydro¬ 
gen,  for  a  given  value  of  power  demand.  Based  on  the  formulation  and  solution  of  the  NLP  problem,  a 
look-up  table  is  developed,  containing  the  optimal  values  of  the  system  variables  for  any  possible  value 
of  power  demand.  In  the  last  step,  a  Model  Predictive  Control  (MPC)  methodology  is  designed,  for  the 
optimal  control  of  the  system  response  to  successive  sep-point  changes  of  power  demand.  The  efficiency 
of  the  produced  MPC  system  is  illustrated  through  a  number  of  simulations,  which  show  that  a  successful 
dynamic  closed-loop  behaviour  can  be  achieved,  while  at  the  same  time  the  consumption  of  hydrogen  is 
minimized. 

©  2009  Elsevier  B.V.  All  rights  reserved. 


1.  Introduction 

A  fuel-cell  system  is  a  device  where  chemical  energy  from  a 
fuel,  such  as  hydrogen,  is  electrochemically  converted  to  electri¬ 
cal  and  thermal  energy,  without  combustion.  Fuel  cells  operate 
at  low  noise  levels  and  do  not  produce  pollutant  emissions. 
Fuel  cells  offer  various  economic  and  environmental  advantages 
over  internal-combustion  engines  and  batteries.  Due  to  those 
advantages,  fuel-cell  technology  has  come  to  the  foreground  for 
many  mobile  and  terrestrial  applications.  In  particular,  the  Pro¬ 
ton  Exchange  Membrane  Fuel-Cell  (PEMFC)  technology  offers  many 
advantages  over  other  types  of  fuel  cells  (FCs);  hence,  this  device 


Abbreviations:  3D,  Three-Dimensional;  CPU,  Central  Processor  Unit;  DC,  Direct 
Current;  FIR,  finite  impulse  response;  FC,  fuel  cell;  LOO,  Leave  One-Out;  LQG, 
Linear  Quadratic  Gaussian;  LQR,  Linear  Quadratic  Regulator;  MIMO,  multi-input- 
multi-output;  MPC,  Model  Predictive  Control;  NLP,  Non-Linear  Programming; 
NNM,  neural  network  model;  PEMFC,  Proton  Exchange  Membrane  Fuel-Cell;  PID, 
proportional-integral-derivative;  PRESS,  Prediction  Error  Sum  of  Squares;  RBF, 
Radial  Basis  Function;  RMSE,  Root  Mean  Squared  Error;  SOFC,  Solid  Oxide  Fuel  Cell; 
SSE,  Sum  of  Squared  Errors  between  the  observations  and  the  predicted  values;  SSY, 
Sum  of  Squared  Deviations  between  the  observations  and  their  mean. 
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has  the  potential  to  become  a  primary  power  source  for  power  in 
the  coming  era. 

However,  the  FC  technology  needs  to  tackle  a  variety  of  signifi¬ 
cant  problems,  in  order  to  become  more  competitive  and  efficient 
compared  to  other  energy  technologies.  In  particular,  fuel-cell  sys¬ 
tems  pose  many  challenging  problems  as  far  as  control  and  system 
integration  are  concerned,  due  to  complicated  system  configura¬ 
tions  and  interactive  subsystems,  which  lead  to  intricate  transient 
dynamics  and  difficult  control  tasks. 

In  the  literature,  a  number  of  control  strategies  and  method¬ 
ologies  for  fuel-cell  systems  have  been  proposed,  that  range 
from  simple  proportional-integral-derivative  (PID)  controllers  to 
advanced  control  strategies  such  as  fuzzy  controllers  and  Model 
Predictive  Control  (MPC)  methodologies. 

Methekar  et  al.  [1]  considered  a  multi-input-multi-output 
(MIMO)  system  with  inputs  of  hydrogen  and  coolant  and  outputs 
of  power  density  and  temperature,  and  proposed  two  PID  control 
strategies.  Woo  and  Benziger  [2]  demonstrated  that  the  power  out¬ 
put  from  a  PEMFC  can  be  directly  regulated  by  limiting  the  hydrogen 
feed  to  the  FC.  Regulation  was  accomplished  by  varying  the  internal 
resistance  of  the  membrane-electrode  assembly  in  a  self-draining 
FC  with  the  effluents  connected  to  water  reservoirs. 

Pukrushpan  et  al.  [3]  proposed  both  a  feedforward  controller 
and  a  Linear  Quadratic  Regulator  (LQR)  state  feedback  controller,  for 
the  online  control  of  a  PEMFC  system.  Guidong  et  al.  [4]  presented 
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Nomenclature 

ch  control  horizon 

C  number  of  past  inputs  in  the  FIR  model 

d  correction  vector 

/  Radial  Basis  Function 

g\ ,  gi  function  symbols 

h[  ith  coefficient  in  FIR  model 

Hj  matrix  containing  the  ith  coefficients  of  a  MIMO  FIR 

model 

1st  stack  current  (A) 

I*t  optimal  value  of  stack  current  (A) 

k  number  of  independent  variables 

I  number  of  hidden  nodes  in  the  neural  network 

topology 

M  number  of  output  variables 

n  number  of  the  available  data 

N  number  of  input  variables 

ph  prediction  horizon 

PD  power  demand  (W) 

Pnet  fuel-cell  stack  net  power 

R  weight  matrix 

R2  coefficient  of  determination  (%) 

R2v  coefficient  of  determination  by  using  the  cross- 

validation  technique  (%) 

Tt  sampling  time  (s) 

u  input  vector 

Au  vector  of  control  moves 

Ucm  compressor  voltage  (V) 

Ust  stack  voltage  (V) 

U*m  optimal  value  of  compressor  voltage  (V) 

Wj  the  weight  corresponding  to  the  response  of  the  jth 
node  in  the  neural  network  topology 
W  weight  matrix 

x  input  vector  in  the  neural  network  model 

Xj  the  centre  of  the  jth  node  in  the  neural  network 

topology 

y  output  vector 

y  prediction  of  output  vector 

ym  i  true  value  of  the  mth  output  variable  for  the  ith 
observation 

ym  i  predicted  value  of  the  mth  output  variable  for  the 
ith  observation 

ym  mean  of  all  values  of  the  mth  output  variable  over 
the  available  data  set 

pLoo  rnth  output  variable  prediction  for  observation  i  of 
the  model  that  is  trained  using  all  the  available  data, 
except  from  observation  i 
y sp  set-point  vector 

Zj  The  response  of  the  jth  node  in  the  neural  network 

topology 

Greek  symbols 

Xq2  oxygen  excess  ratio 

A*  optimal  value  of  oxygen  excess  ratio 


a  power  management  system,  by  employing  a  DC-DC  module,  a 
battery,  a  controller  and  some  other  circuits,  as  well  as  an  LQR 
algorithm,  which  provided  load  with  an  approximately  constant 
voltage,  but  also  stabilized  stack  current  by  regulating  charging  or 
discharging  current  of  the  battery,  according  to  the  change  of  load. 
Rodatz  et  al.  [5]  designed  a  Linear  Quadratic  Gaussian  (LQG)  con¬ 
troller  to  decouple  the  pressure  trace  from  the  mass  flow  trace.  Di 


Domenico  et  al.  [6]  extended  this  idea,  by  developing  a  multivari¬ 
able  LQG  controller  that  is  able  to  tune  the  excess  air  ratio  while 
tracking  the  optimal  pressurization  to  maximize  system  efficiency 
for  transient  loads. 

Danzer  et  al.  [7]  proposed  a  multivariable  flatness-based  control 
technique,  for  controlling  cathode  pressure  and  oxygen  excess  ratio 
of  a  fuel  cell,  using  the  mass  flow  controller  and  the  outlet  throttle 
as  actuators.  Zenith  and  Skogestad  [8]  presented  a  method  to  con¬ 
trol  the  output  of  a  buck-boost  converter  connected  to  a  fuel  cell, 
by  switching  the  converter  using  a  few  logical  rules.  The  resulting 
control  loop  was  then  inserted  in  a  cascade-control  framework  to 
control  the  armature  current  in  a  DC  motor,  by  manipulating  its 
input  voltage. 

Wang  et  al.  [9]  applied  multivariable  robust  control  strate¬ 
gies  to  a  PEMFC  system.  Initially,  they  modeled  the  PEMFC  as  a 
two-input-two-output  system,  where  the  inputs  are  air  and  hydro¬ 
gen  flow  rates  and  the  outputs  are  cell  voltage  and  current.  By 
fixing  the  output  resistance,  they  aimed  to  control  the  cell  volt¬ 
age  output  by  regulating  the  air  and  hydrogen  flow  rates.  An 
multivariable  robust  controller  was  designed  to  provide  robust 
performance  and  reduce  the  hydrogen  consumption  of  this  sys¬ 
tem. 

Mo  Zhijun  et  al.  [10]  applied  a  fuzzy  control  technique  to  a 
PEMFC,  in  order  to  maintain  the  output  voltage  through  the  gas 
pressure  control,  for  any  load  applied  to  the  PEMFC  output  ter¬ 
minals.  Zhan  Yuedong  et  al.  [11]  designed  a  fuzzy-PI  controller  to 
control  the  hydrogen  and  air/oxygen  mass  flows,  but  also  auxiliary 
variables  such  as  the  temperature,  pressure,  humidity  of  the  mem¬ 
brane,  and  proportion  of  stoichiometry.  Sakhare  et  al.  [12]  designed 
a  power  condition  unit  for  a  Solid  Oxide  Fuel  Cell  (SOFC)  by  incor¬ 
porating  a  fuzzy  logic  control  strategy,  which  is  also  suitable  for 
any  type  of  fuel  cells  with  slight  alterations,  whereas  Schumacher 
et  al.  [13]  proposed  a  scheme  for  the  online  control  of  multiple 
cell  units,  using  the  Takagi  Sugeno  Kang  (TSK)  type  fuzzy  logic 
approach. 

Arce  et  al.  [14]  developed  an  explicit  predictive  control  strat¬ 
egy  for  a  stand-alone  PEMFC.  Golber  and  Lewin  [15]  presented  an 
alternative  model-based  controller  for  the  regulation  of  a  PEMFC. 
The  model  accounts  for  spatial  dependencies  of  voltage,  current, 
material  flows,  and  temperatures  in  the  fuel-cell  channel.  Analysis 
of  the  process  model  showed  that  the  effective  gain  of  the  process 
undergoes  a  sign  change  in  the  normal  operating  range  of  the  fuel 
cell,  indicating  that  it  cannot  be  stabilized  using  a  linear  controller 
with  integral  action  and,  consequently,  they  developed  a  non-linear 
model  predictive  controller.  Wang  et  al.  [16]  applied  a  data-driven 
predictive  control  approach  in  order  to  solve  the  control  problem 
of  a  SOFC  system,  since  it  was  previously  shown  that  the  control 
of  SOFC  is  challenging,  due  to  the  slow  response  and  strict  oper¬ 
ating  constraints.  An  MPC  strategy  was  developed  by  Vahidi  et  al. 
[17]  for  a  hybrid  PEMFC  system  with  ultracapacitors  as  an  auxiliary 
source  of  power.  Jurado  [18]  developed  an  MPC  method  for  SOFCs, 
by  integrating  in  it  fuzzy  models.  The  MPC  concept  was  also  used 
by  Vahidi  et  al.  [19]  to  avoid  fuel-cell  oxygen  starvation,  prevent  air 
compressor  surge  and  choke,  and  match  at  the  same  time  an  arbi¬ 
trary  level  of  current  demand.  The  control  scheme  was  designed 
on  a  hybrid  configuration,  in  which  a  bank  of  ultracapacitors  sup¬ 
plements  the  polymer  electrolyte  membrane  fuel  cell  during  fast 
current  transients. 

In  this  paper,  a  new  methodology  is  presented  for  the  opti¬ 
mal  control  of  a  PEMFC  system.  A  simulation  model  is  used  first 
to  generate  a  database  that  contains  the  values  of  the  key  system 
variables:  the  FC  current,  the  compressor  motor  voltage,  the  FC  net 
power  and  the  oxygen  excess  ratio.  An  advanced  neural  network 
model  (NNM)  methodology  is  then  applied  to  obtain  a  non-linear 
model  between  the  input  and  the  output  variables.  The  database 
is  also  used  for  validating  the  accuracy  of  the  produced  models. 
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Fig.  1.  Flowchart  of  the  proposed  PEMFC  modeling,  optimization  and  control  pro¬ 
cedure. 


Table  1 

Structure  of  the  database  used  for  the  generation  and  validation  of  the  meta-models. 


No.  of  data 

Input  variables 

Output  variables 

ht  (A) 

Ucm  (V) 

Pnet  (W) 

^o2 

1 

60 

50 

14,580 

1.312 

2 

100 

65 

21,640 

1.176 

3 

150 

90 

29,750 

1.235 

75 

60 

17,600 

1.369 

524 

P-f  (1st,  Vcm) 


200 


150 


100 


1st  A 


250 


x  10 


150 


Vcm  (V) 


Fig.  2.  FC  stack  net  power  as  a  function  of  stack  current  and  compressor  voltage. 


2.2.  Formulation  of  database 


A  Non-Linear  Programming  (NLP)  optimization  problem  is  formu¬ 
lated  next,  in  order  to  obtain  the  optimal  values  of  the  decision 
variables.  Eventually,  a  methodology  based  on  the  principles  of  MPC 
is  designed,  for  the  optimal  control  of  the  system  response  to  suc¬ 
cessive  sep-point  changes  of  power  demand.  The  optimal  oxygen 
excess  ratio  value  (determined  by  the  solution  of  the  optimization 
problem  as  described  above)  is  used  as  an  additional  controlled 
variable.  The  efficiency  of  the  produced  MPC  system  is  illustrated 
on  a  number  of  simulations,  which  show  that  a  successful  dynamic 
closed-loop  behaviour  can  be  achieved,  while  at  the  same  time  the 
consumption  of  hydrogen  is  minimized. 

The  proposed  computational  tool  is  summarized  graphically  in 
Fig.  1. 


The  available  first-principle  model  was  used  to  provide  the 
steady-state  output  values  corresponding  to  a  number  of  input 
values  that  cover  the  full  operational  range  of  the  system.  In  par¬ 
ticular,  values  in  the  range  (20A-350A)  were  assigned  to  stack 
current,  whereas  compressor  voltage  was  set  at  values  in  the  range 
(20  V-250  V).  CPU  time  for  each  run  was  about  1  s  on  an  Intel  Core 
2.2  GHz  processor.  The  obtained  results  formed  a  database,  which 
consists  of  524  input-output  examples.  Table  1  presents  the  struc¬ 
ture  of  the  database  and  includes  some  indicative  values  of  the 
input-output  variables. 

The  contents  of  the  database  are  presented  graphically  in 
Figs.  2  and  3.  Fig.  2  presents  in  3D  view  the  net  power  as  a  func¬ 
tion  of  stack  current  and  compressor  voltage,  whereas  Fig.  3  shows 


2.  PEMFC  modeling 

2.1.  Dynamic  first-principle  PEMFC  model 

The  four-step  methodology  that  is  presented  in  this  paper  is  built 
upon  a  detailed  first-principle  dynamic  model  of  a  PEMFC  system, 
which  has  been  proposed  by  Pukrushpan  et  al.  [3].  The  model  con¬ 
tains  both  the  main  part  (cell  units)  and  auxiliary  parts  (compressor, 
cooling  system,  humidification  system  as  well  as  the  supply  mani¬ 
folds)  of  the  PEMFC  system.  The  overall  inputs  to  the  model  are  the 
stack  current,  Ist,  and  the  compressor  voltage,  Ucm,  whereas  model 
outputs  are  the  FC  stack  net  power,  Pnetf  and  the  oxygen  excess  ratio, 
Xq2,  (input  oxygen  mass  rate/oxygen  reacted).  More  details  about 
the  development  of  the  model,  the  complete  set  of  differential  and 
algebraic  equations  and  the  values  of  the  model  parameters  are 
given  in  the  above  reference. 


X02=f  (Ist,  Vcm) 


Fig.  3.  Oxygen  excess  ratio  as  a  function  of  stack  current  and  compressor  voltage. 
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INPUT  1  INPUT  2  INPUT  N 

Fig.  4.  Standard  RBF  neural  network  topology. 


the  oxygen  excess  ratio,  as  a  function  of  the  same  input  variables.  It 
should  be  noted  that  in  case  we  have  access  to  experimental  data 
obtained  from  a  true  fuel-cell  system,  these  data  can  be  stored  in  the 
database,  instead  of  the  first-principle  model  simulation  results. 

3.  Development  of  meta-models 

The  database  that  was  generated  as  described  in  the  previous 
section  was  utilized  for  developing  correlation  equations  between 
the  output  variables  Pnet,  and  Xq2  and  the  input  variables  Ist  and  UCm 
of  the  system.  In  particular,  the  Radial  Basis  Function  (RBF)  neural 
network  architecture  was  adopted  for  producing  the  meta-models, 
due  to  its  simple  topology  and  the  fast  and  robust  algorithms  that 
are  available  in  the  literature  for  training  this  type  of  networks.  The 
topology  of  the  RBF  network  is  presented  in  Fig.  4  and  consists  of 
three  layers:  the  input  layer,  the  hidden  layer  and  the  output  layer. 

The  input  vector  to  the  network  u  contains  the  input  variables, 
i.e.u  =  [«!,  u2]t  =  [1st,  UCm]T,  while  the  neural  network  output  vec¬ 
tor  provides  the  estimated  values  for  the  net  power  and  oxygen 

excess  ratio  y  =  [y\,y2^  =  [Pnet,^ o2]  -  which  are  calculated  as 
weighted  summations  of  the  hidden  node  responses: 

L 

ym  -  y^wm j  ■  zj(x)  (l) 

j= 1 

where 

Zj(u)=/(l|u-Xj|||)  (2) 

In  the  above  equations  Zj  is  the  response  of  the  jth  node,  /  is 
the  Radial  Basis  Function,  Xj  is  the  centre  of  the  jth  node,  L  is  the 
total  number  of  hidden  nodes  and  wmj-  is  the  weight  connecting 
the  response  of  the  jth  node  to  the  mth  network  output. 

An  RBF  training  procedure  aims  at  the  determination  of  the 
number  of  nodes  in  the  hidden  layer,  the  hidden  node  centres  and 
the  output  weights,  in  order  to  minimize  the  deviation  between  the 
predicted  and  the  measured  values  of  the  output  variables  over  the 
set  of  the  available  data  base.  The  training  method  used  in  this  work 
is  based  on  the  fuzzy  partition  of  the  input  space,  which  is  produced 
by  defining  a  number  of  triangular  fuzzy  sets  in  the  domain  of  each 
input  variable  [20].  The  centres  of  these  fuzzy  sets  form  a  multidi¬ 
mensional  grid  on  the  input  space.  A  rigorous  selection  algorithm 
chooses  the  most  appropriate  vertices  on  the  grid,  which  are  then 
used  as  the  hidden  node  centres  in  the  resulting  RBF  network  model. 
The  idea  behind  the  selection  algorithm  is  to  place  the  centres  in 


the  multidimensional  input  space,  so  that  the  distance  between 
any  two  centre  locations  is  guaranteed  to  be  greater  than  a  lower 
limit,  which  is  defined  by  the  length  of  the  edges  on  the  grid.  At  the 
same  time,  the  algorithm  ensures  that  for  any  input  example  in  the 
training  set  there  is  at  least  one  selected  hidden  node  that  is  close 
enough,  according  to  an  appropriately  defined  distance  criterion. 
The  so  called  “fuzzy-means”  training  method  does  not  need  the 
number  of  centres  to  be  fixed  before  the  execution  of  the  method. 
Due  to  the  fact  that  it  is  a  one-pass  algorithm,  it  is  extremely  fast, 
even  in  the  case  of  a  large  database  of  input-output  training  data. 
One  additional  advantage  is  that  the  training  algorithm  used  needs 
only  one  tuning  parameter,  namely  the  number  of  fuzzy  sets  that 
are  utilized  to  partition  each  input  dimension. 


3.1  Correlation  equations 


The  training  procedure  was  used  several  times  by  altering  each 
time  the  fuzzy  partition  of  the  input  space  (number  of  fuzzy  sets 
defined  in  each  input  dimension),  which  is  in  fact  the  only  design 
parameter  that  must  be  defined  by  the  user,  when  utilizing  the 
fuzzy-means  algorithm.  The  average  time  needed  to  train  an  RBF 
network  using  an  Intel  Core  2.2  GFIz  processor  was  18  s.  The  results 
were  gradually  improved  up  to  the  point  where  16  fuzzy  sets  were 
used  to  partition  the  domain  of  each  input  variable,  which  cor¬ 
responds  to  57  hidden  nodes.  A  further  increase  results  to  the 
overtraining  phenomenon,  where  the  performance  of  the  produced 
model  is  not  improved,  although  the  model  increases  in  size.  The 
following  statistics  were  used  to  measure  the  accuracy  of  the  pro¬ 
duced  models: 


Root  Mean  Squared  Error 


RMSEm  = 


SSErr 


n-k-  1 


TVi  ( ym.i  ym,i ) 

n-k  -  1 


(3) 


Coefficient  of  Determination 

p2  _  -i  SSEm  ^j_i  {ym.i  ~  9mj) 

m~  “  £"=1(3v.--ym)2 


(4) 


F-statisticm 


((1-R2  )/(n -/<-!)) 


(5) 


where  SSEm  measures  the  Sum  of  Squared  Errors  between  the  obser¬ 
vations  and  the  predicted  values  over  the  set  of  the  available  data  for 
the  mth  output  variable,  while  SSYm  measures  the  respective  sum¬ 
mation  of  squared  deviations  between  the  observations  and  their 
mean.  In  the  above  equations  ym  i  is  the  true  value  of  the  mth  output 
variable  for  the  ith  observation,  ym  i  is  the  associated  model  predic¬ 
tion,  and  ym  is  the  mean  of  all  values  of  the  mth  output  variable  in 
the  available  data  set.  Finally,  k  is  the  number  of  the  independent 
variables,  and  n  the  number  of  the  available  input-output  data. 

In  order  to  test  the  reliability  of  the  modeling  methodology, 
the  cross-validation  method  was  used.  Following  this  technique, 
a  number  of  modified  data  sets  is  created  by  deleting  in  each  case 
one  or  a  small  group  (i.e.  leave  some  out)  of  objects.  For  each  data 
set,  an  input-output  model  is  developed,  by  applying  the  model¬ 
ing  technique.  Each  model  is  evaluated,  by  measuring  its  accuracy 
in  predicting  the  responses  of  the  remaining  data  (the  ones  that 
have  not  been  utilized  in  the  development  of  the  model).  In  par¬ 
ticular,  the  Leave  One-Out  (LOO)  cross-validation  procedure  was 
used  in  this  study,  which  produces  a  number  of  models,  by  deleting 
each  time  one  object  from  the  training  set.  Obviously,  the  number 
of  models  produced  by  the  LOO  procedure  is  equal  to  the  number 
of  available  examples  n.  Prediction  Error  Sum  of  Squares  {PRESS)  is 
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Table  2 

Statistical  indices  corresponding  to  the  produced  RBF  neural  network  model  for  both 
output  variables. 


Statistical  indices  relating  to  RBF  neural 
network  model  with  57  hidden  nodes 

Pnet 

A-o2 

RMSE 

52.9235 

0.0426 

R2 

0.9999 

0.9997 

F-Statistic 

2,604,739.5 

867,812.4 

R2CV 

0.9998 

0.9996 

SpRESS 

98.2385 

0.0513 

a  standard  index  to  measure  the  accuracy  of  a  modeling  method 
using  the  cross-validation  technique.  Based  on  the  PRESS  and  SSY 
statistics,  the  R2V  and  SPREss  values  can  be  easily  calculated  for  each 
output  variable.  The  formulae  used  to  calculate  all  the  aforemen¬ 
tioned  statistics  are  presented  below: 


p2 

KCV,m 


=  1  - 


PRESSm 

~SSY^~ 


E"=i(yn,.-p2 
ELi  (3Vi  -  yf 


(6) 


Fig.  6.  Database  oxygen  excess  ratio  values  against  neural  network  predictions  using 
the  cross-validation  procedure. 


$  PRESS,  m  = 


PRESSm 


n-k-  1 


(7) 


where  y is  the  mth  output  variable  prediction  for  observation 
i  given  by  the  model  that  has  been  trained  using  all  the  available 
data,  except  from  observation  i. 

The  statistical  indices  corresponding  to  the  produced  RBF  neural 
network  model  for  both  output  variables  are  summarized  in  Table  2. 

Figs.  5  and  6  depict  the  true  FC  stack  net  power  and  oxygen  excess 
ratio  values,  against  the  ones  produced  by  the  RBF  neural  network 
model,  using  the  LOO  cross-validation  procedure.  Good  predictive 
abilities  of  both  models  are  illustrated  by  noticing  that  in  both  cases 
all  points  lie  very  close  to  the  main  diagonal  of  the  graph. 

The  RBF  neural  network  modeling  technique  was  further  val¬ 
idated  by  developing  models  based  on  only  75%  of  the  available 
data  and  evaluating  them  on  the  validation  set,  i.e.  the  rest  of  the 
database  that  was  not  used  for  the  derivation  of  the  model.  In  order 
to  show  that  the  success  of  one  particular  model  is  not  due  to  a 
chance  correlation,  100  random  partitions  of  the  data  into  train¬ 
ing  and  validation  sets  (75  and  25%  of  the  data,  respectively)  were 
performed  to  generate  100  different  RBF  models.  The  RMSE  and  R2 
statistics  were  calculated  for  each  model  using  only  the  validation 
examples.  Table  3  shows  the  min,  max,  mean  and  standard  devia- 


Table  3 

Min,  max,  mean  and  standard  deviation  corresponding  to  RMSE  and  R2  values 
obtained  after  performing  100  random  partitions  of  the  available  data  into  training 
and  validations  sets. 


tion  of  the  RMSE  and  R2  values  obtained  after  performing  the  100 
runs. 

Fig.  7  presents  the  net  power  database  values  and  the  predic¬ 
tions  of  the  neural  network  associated  with  the  following  statistics 
as  far  as  the  validation  examples  are  concerned:  RMSE=  108.49 
and  R2  =  0.9999.  Accordingly,  Fig.  8  plots  the  oxygen  excess  ratio 
database  values  and  the  predictions  of  the  neural  network  associ¬ 
ated  with  the  following  statistics  as  far  as  the  validation  examples 
are  concerned:  RMSE  =  0.0458  and  R2  =  0.9996.  In  both  cases  the  dif¬ 
ference  between  the  two  lines  is  almost  invisible,  thus  illustrating 
the  predictive  ability  of  the  produced  models. 

It  must  be  mentioned  here  that  before  applying  the  RBF  method¬ 
ology,  several  models  were  developed  using  the  linear  regression 


Fig.  5.  Database  net  power  values  against  neural  network  predictions  using  the 
cross-validation  procedure. 


Fig.  7.  Database  net  power  values  against  neural  network  predictions  for  the  vali¬ 
dation  data  set. 
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Validation  data  set 


Fig.  8.  Database  oxygen  excess  ratio  values  against  neural  network  predictions  for 
the  validation  data  set. 


approach.  The  models  were  tested  using  the  same  validation  rou¬ 
tines  that  were  used  to  test  the  RBF  models.  Best-fitting  equations 
concerning  the  FC  stack  net  power  and  the  oxygen  excess  ratio  are 
given  below: 


Pnet  =  159.15/st  -  1198.59^+7051.57  (8) 

ht 

X0.  =  0.0053^  +  1.3187^+0.1067  (9) 

1st  1st 

The  aforementioned  equations  are  associated  with  the  statistics 
shown  in  Table  4. 

A  comparison  of  the  results  included  in  Tables  3  and  4  clearly 
indicates  that  the  RBF  modeling  methodology  outperforms  the  lin¬ 
ear  regression  approach,  especially  when  prediction  of  net  power 
is  concerned.  The  RBF  model  will  thus  be  used  in  the  formulation 
of  the  optimization  problem,  which  is  described  in  the  next  section 
of  the  paper. 


4.  Optimization  of  the  PEMFC  performance 

Performance  of  fuel-cell  systems  can  be  optimized  by  selecting 
appropriate  stack  current  and  compressor  voltage  input  values  [3]. 
These  are  not  fixed,  but  depend  on  the  power  that  is  produced  by 
the  system.  In  this  part  of  the  paper,  the  meta-model  produced  in 
the  previous  section  is  utilized  to  compute  the  optimal  values  of 
the  input  variables  for  different  levels  of  power  demand.  The  opti¬ 
mization  problem  is  formulated  so  that  it  takes  into  account  the 
various  constraints  and  limitations  of  the  system.  In  Section  4.1  the 
formulation  of  the  optimization  problem  is  presented  in  detail. 


Table  4 

Statistical  indices  corresponding  to  the  best  models  obtained  using  the  linear  regres¬ 
sion  technique. 


Statistical  indices  related  to  the  best  models 
obtained  using  the  linear  regression  technique 

Pnet 

A.o2 

RMSE 

2530.6 

0.1282 

R2 

0.9505 

0.9975 

F-Statistic 

4992.5 

103,740 

R2CV 

0.9496 

0.9975 

SpRESS 

2550 

0.1299 

1 

^ Power  demand 


Decision 

Variables 

1 

^Stack  current 
^Compressor  voltage 
^Oxygen  excess  ratio 


Objective  Function  to 
minimize 


1 


^Stack  current 


Formulation  and  solution 
of  a  Non-Linear  Programming 
Problem 


Fig.  9.  Formulation  of  the  Non-Linear  Programming  Problem. 


4.1.  Formulation  of  the  optimization  problem 

The  optimization  problem  is  completely  defined  by  the  con¬ 
straints  that  should  be  satisfied  by  the  solution  of  the  problem  and 
the  performance  criterion,  which  is  minimization  of  the  FC  stack 
current  for  a  given  load. 

4.1.1.  Constraints 

i)  The  FC  stack  net  power  has  been  expressed  as  a  function  of  Jst, 
Ucm  using  the  RBF  modeling  methodology: 

Pnet  =  gliht,  CfCm)  (10) 

ii)  Accordingly,  an  RBF  model  has  been  developed  that  expresses 
oxygen  excess  ratio  as  a  function  of  Jst,  Ucm‘- 

^02  =g2{Ist,Ucm)  (11) 

iii)  FC  stack  net  power  should  match  the  power  demand  by  the 
consumer  (set-point): 


Pnet  =  Pd  (12) 

A  number  of  additional  constraints  are  used  in  the  optimiza¬ 
tion  problem  to  express  physical  limitations  of  the  system  and/or 
bound  the  system  variables  between  desired  upper  and  lower  lim¬ 
its.  In  particular,  the  following  upper  and  lower  limits  are  imposed: 
net  power  (5kW-45kW),  oxygen  excess  ratio  (>1),  stack  current 
(20  A-350  A),  compressor  voltage  (20V-250V). 

4.1.2.  Objective  function 

Selection  of  the  objective  function  is  a  critical  issue  in  the  formu¬ 
lation  of  the  optimization  problem.  Minimization  of  stack  current 
1st  was  selected  as  the  performance  criterion,  because  due  to  Fara¬ 
day’s  First  Law,  stack  current  is  proportional  to  the  consumption  of 
the  costly  fuel.  Thus,  by  minimizing  stack  current,  consumption  of 
hydrogen  is  also  reduced.  Consequently,  the  optimization  problem 
is  formulated  as  follows: 

min  1st  (13) 

ht, Ucm, Pnet, ^  02 

Subject  to  the  constraints  (10)-(12)  and  the  upper  and  lower 
bounds  posed  on  the  variables  of  the  problem.  The  optimization 
problem  is  presented  graphically  in  Fig.  9. 

The  previously  described  optimization  problem  is  an  NLP  due 
to  the  non-linear  neural  network  constraints  (10)  and  (11).  The 
global  minimum  of  the  problem  was  obtained  using  the  GAMS 
modeling  language  [21  ].  The  solution  of  the  optimization  problem 
contains  the  optimal  values  of  all  decision  variables  (optimal  values 
for  the  manipulated  variables  stack  current  and  compressor  volt¬ 
age,  I*t,  U*m,  as  well  as  the  optimal  values  for  oxygen  excess  ratio, 


264 


J.  Hasikos  et  al.  /  Journal  of  Power  Sources  193  (2009)  258-268 


Fig.  10.  Optimal  A,*  ^  -values  as  a  function  of  power  demand. 


4.2.  Solution  of  the  optimization  problem— results 

The  previously  described  optimization  tool  was  used  to  obtain 
the  optimal  values  of  stack  current,  compressor  voltage  and  oxygen 
excess  ratio,  corresponding  to  several  levels  of  power  demand,  that 
cover  the  range  of  operation  (5kW-45kW).  The  results  are  pre¬ 
sented  graphically  in  Figs.  10  and  11.  In  particular,  Fig.  10  depicts 
the  optimal  values  of  the  ratio  A*  as  a  function  of  power  demand, 
whereas  Fig.  11  illustrates  the  optimal  values  of  the  input  variables 
I*t,  U*m.  These  figures  will  prove  useful  when  designing  the  process 
control  strategies  in  the  next  section  of  this  paper. 


5.  Control  of  the  PEMFC 

The  final  step  in  the  development  of  the  proposed  framework 
is  to  design  the  optimal  control  strategy  for  the  FC  system.  Control 
of  fuel  cells  is  a  challenging  problem,  due  to  the  multiple  and  con¬ 
tradicting  criteria  that  define  the  performance  of  the  system.  The 
first  priority  in  a  fuel-cell  control  system  is  to  respond  both  quickly 
and  smoothly  to  power  demand  changes.  This  is  particularly  impor¬ 
tant  for  applications  involving  frequent  changes  of  power  demand, 
such  as  in  automobile  systems.  However,  changes  in  power  demand 
cause  sudden  changes  in  the  oxygen  excess  ratio  values  that  are 
often  responsible  for  the  oxygen  depletion  phenomenon,  which  in 
turn  can  lead  to  membrane  damage  and  to  an  abrupt  FC  system 
shut-down. 


The  MPC  methodology  was  chosen  to  design  a  control  strategy 
which  improves  the  dynamic  behaviour  of  the  system,  while  at  the 
same  time  the  fuel  gas  consumption  is  reduced.  MPC  has  advanced 
to  a  popular  control  methodology  for  industrial  and  process  appli¬ 
cations  mainly  due  to  the  inherent  ability  of  the  method  to  handle 
efficiently  constraints  and  uncertainties  in  multivariable  dynami¬ 
cal  systems.  Specifically,  the  dynamic  matrix  control  (DMC)  method 
was  adopted  (Cutler  and  Ramaker  [22]),  which  requires  simple  step- 
(or  pulse-)  discrete-time  response  models,  that  can  easily  represent 
time  delays  and  complex  dynamics.  Using  the  DMC  methodology,  at 
each  sampling  instance,  the  sequence  of  the  manipulated  variable 
values  over  a  future  control  horizon  is  computed,  by  solving  online 
an  optimization  problem  which  minimizes  the  offset  (deviation  of 
the  controlled  variables  from  the  desired  set-points)  over  a  future 
prediction  horizon,  and  the  control  energy  over  a  future  control 
horizon.  The  disturbance  during  the  prediction  horizon  is  assumed 
to  be  constant  and  equal  to  the  difference  between  the  process 
output  and  the  model  prediction  at  the  current  time  instance. 


5.1.  Design  of  a  DMC  methodology 


5.1.1.  Development  of  finite  impulse  response  (FIR)  models 

An  FIR  or  convolution  model  is  a  discrete  time  filter  that  corre¬ 
lates  dynamically  an  output  variable  of  the  system  with  an  input 
variable,  as  follows: 

c 

y(t)  =  ^hi.u(t-i)  (14) 

i= 1 


where  C  is  the  number  of  past  inputs  used  in  the  FIR  model.  The 
coefficients  in  the  above  model  can  be  easily  determined  from  the 
dynamic  response  of  the  output  variable,  after  the  introduction  of 
a  step  or  pulse  change  on  the  input  variable,  when  the  system  is 
initially  at  steady  state  (Seborg  et  al.  [23  ]).  In  case  of  a  MIMO  system, 
separate  FIR  models  are  developed  for  each  pair  of  input-output 
variables.  The  FIR  models  can  however  be  aggregated  in  a  matrix 
form  as  follows: 


y  (/<)  = 


y  i(/<) 
y2(fc) 


=  Y.Hi ' u(t 

i=i 


0 


=  £ 


h]'2 

fc2,2 


hM,2 


'ui(fc-i)' 

u2(k  -  i) 

_uN{k-i)_ 

(15) 


where  h™,n  is  the  ith  coefficient  in  the  FIR  model  correlating  the 
mth  output  variable  with  the  nth  input  variable. 

The  fuel-cell  system  of  the  present  study  is  a  2  x  2  MIMO  sys¬ 
tem.  The  FIR  models  were  produced  by  applying  step  changes  on 
the  manipulated  variables  ( Jst,  Ucm )  when  the  system  was  initially  at 
steady  state,  corresponding  to  Jst=100A,  Ucm  =99  V  (the  respective 
steady-state  output  values  are  Pnet  =  23,255  W,  Aq2  =  2).  The  sam¬ 
pling  time  was  set  to  0.075  s,  so  that  the  FIR  model  captures  the 
dynamic  characteristics  of  the  system,  but  at  the  same  time  uses  a 
reasonable  number  of  past  input  values  (C= 30). 


5.1.2.  The  DMC  algorithm 

In  the  DMC  methodology,  at  time  instance  k  a  set  of  future 
manipulated  variables  (sequence  of  manipulated  variable  values 
over  the  control  horizon)  is  selected,  in  order  to  minimize  an  appro¬ 
priate  objective  function.  The  objective  function  includes  both  the 
deviations  of  the  predicted  controlled  variables  from  their  set- 
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Fig.  12.  DMC  basic  concept. 


points  over  the  future  prediction  horizon  (ph),  and  the  control  effort 
over  the  control  horizon  (ch).  Then,  only  the  first  control  move 
is  implemented  and  the  optimization  problem  is  formulated  and 
solved  again  at  the  next  time  step  k  + 1.  Fig.  12  presents  a  schematic 
representation  of  the  DMC  philosophy. 

In  the  sequel,  the  optimization  problem  that  is  solved  at  each 
time  instance  is  presented  in  detail.  It  consists  of  the  objective  func¬ 
tion,  (Eq.  (16))  which  describes  both  the  aforementioned  control 
targets  and  a  number  of  constraints,  (Eqs.  (17)— (21 )): 

Objective  function: 

min 

Au(k),  Au(k+1 Au(k+ch) 


ph 

y~J|W(y(fe  +  i|fc)  —  ysp(fc))ll|  4 

1=1 

subject  to: 

y(k  +  m=  ^H,u(k  +  i-l)j 

d(k)  =  y(k)-y(k\k)  i  =  l,...,ph 
Au (k  +  i)  =  u (k  +  i)  -  u (k  +  i  -  1 ) 
Umin  <  U (k  +  i)  <  Umax  i  —  1  ?  ...  ? 
Au(/<  +  i)  =  0  i  =  ch  +  1 , . . . ,  ph 


ch  \ 

-^J|RAu(k  +  i)ll2  1 
i=0  / 

(16) 

+  d(/<),  z  =  l,...,  ph 

(17) 

(18) 

i  =  1 , . . . ,  ch 

(19) 

ch 

(20) 

(21) 

In  the  above  formulation  of  the  optimization  problem,  y (k  +  i\k ) 
is  the  i-step  (ahead)  prediction  of  the  controlled-variable  vector, 
y Sp{k)  is  a  vector  containing  the  current  set-point  values  of  the  con¬ 
trolled  variables,  Au(/<  +  i)  are  the  future  control  moves,  R,  W  are 
diagonal  square  weight  matrices,  and  d(/<)  is  a  correction  term,  cal¬ 
culated  as  the  current  (at  time  instance  k)  difference  between  the 
true  output  vector  and  the  respective  model  prediction  (Eq.  (18)). 
It  should  be  noted  that  although  the  FIR  model  is  able  to  predict 
the  future  dynamic  behaviour  of  the  fuel-cell  system  (Eq.  (17)),  it 
is  assumed  that  the  future  set-point  changes  and  disturbances  are 
completely  unknown.  The  set-point  vector  is  thus  set  to  its  cur¬ 
rent  value  ysp(/<)  throughout  the  entire  prediction  horizon  and  the 
same  happens  with  the  correction  term  d(/<).  Eq.  (19)  computes  the 
control  move  at  each  time  instance  k  +  i,  as  the  difference  between 
the  input  vector  at  the  same  time  instance  and  the  input  vector  in 
the  previous  time  instance.  Eq.  (20)  contains  hard  constraints  that 


bound  the  manipulated  variables  between  upper  (umax)  and  lower 
(umin)  values  due  to  physical  limitations.  Eq.  (21)  ensures  that  no 
control  moves  are  allowed  after  the  end  of  the  control  horizon. 


5.2.  Implementation  of  the  DMC  controller  on  the  PEMFC  system 

The  DMC  control  strategy  combined  with  the  operational  opti¬ 
mization  results  presented  in  the  previous  section  are  integrated 
in  this  section,  in  order  to  achieve  a  satisfactory  dynamic  perfor¬ 
mance,  while  at  the  same  time  fuel  consumption  is  minimized.  The 
proposed  methodology  is  implemented  as  follows. 

At  each  time  instance  the  set-point  for  the  two  controlled 
variables  are  defined.  The  set-point  for  net  power  demand  is  deter¬ 
mined  by  the  current  energy  needs  of  the  user  of  the  fuel-cell 
system.  The  frequency  of  set-point  changes  in  energy  demand 
depends  on  the  particular  fuel-cell  application.  It  is  most  prob¬ 
able,  that  set-point  changes  will  be  more  frequent  for  a  vehicle 
fuel-cell  system,  compared  to  a  stationary  fuel-cell  application.  As 
mentioned  in  the  formulation  of  the  optimization  problem  in  the 
previous  subsection,  the  control  methodology  does  not  assume  any 
knowledge  about  future  set-point  changes.  The  set-point  for  the 
second  controlled  variable  (oxygen  excess  ratio)  is  then  automati¬ 
cally  determined  by  the  results  presented  in  Section  4.  In  particular, 
the  function  corresponding  to  Fig.  10  is  used  to  obtain  the  oxygen 
excess  ratio  value  that  minimizes  fuel  gas  consumption.  Then,  the 
optimization  problem  defined  by  Eqs.  (16)-(21)  is  solved  to  obtain 
the  optimal  sequence  of  manipulated  variables.  In  the  present  sim¬ 
ulations,  the  prediction  and  control  horizons  were  set  to  ph  =  50 
and  ch  =  40,  respectively.  From  the  optimal  sequence  of  manipulated 
variables,  only  the  first  element  is  actually  applied  to  the  system, 
for  a  time  interval  equal  to  the  sampling  period  that  was  used  to 
obtain  the  FIR  model  coefficients.  At  the  end  of  this  period  the  entire 
procedure  is  repeated. 

The  proposed  methodology  was  tested  on  the  scenario  pre¬ 
sented  in  Table  5  which  involves  a  number  of  step  changes  on  the 
power  demand  set-point.  As  mentioned  in  the  beginning  of  the  sec¬ 
tion,  the  two  performance  criteria  (satisfactory  dynamic  response 
of  both  controlled  variables,  Pnet,  Aq2  )  are  contradictory.  It  will  be 
shown  that  by  tuning  appropriately  the  weight  matrices  in  Eq.  (16), 
special  emphasis  can  be  given  on  one  over  the  other  performance 
criterion.  In  particular: 


Case  I:  The  weight  matrices  are  W  = 


0.003  0 

0  100 


R  = 


0.01  0 
0  1/1.3 

are  shown  in  Figs.  13  and  14. 


and  the  responses  of  the  controlled  variables 


Case  II:  The  weight  matrices  are  W  = 


100  0 
0  10 


R  = 


0.01  0 
0  1/1.3 

are  shown  in  Figs.  15  and  16. 


and  the  responses  of  the  controlled  variables 


Table  5 

A  testing  scenario  for  the  control  methodologies  consisting  of  several  net  power 
set-point  step  tests. 
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Fig.  13.  Response  of  net  power  (DMC,  Case  I). 


Fig.  16.  Response  of  oxygen  excess  ratio  (DMC,  Case  II). 
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Fig.  14.  Response  of  oxygen  excess  ratio  (DMC,  Case  I). 


The  efficiency  of  the  method  is  illustrated  by  observing  that 
the  control  strategy  is  satisfactory  in  both  cases,  because  both 
controlled  variables  respond  quickly  and  the  system  reaches  the 
desired  set-points  with  no-offset.  However,  the  DMC  tuning  of  Case 
I  clearly  favours  the  dynamic  response  of  oxygen  excess  ratio,  which 
is  driven  to  the  desired  set-point  after  some  oscillations  of  very  low 
amplitude,  whereas  1  s  is  needed  for  the  net  power  to  reach  the 
desired  set-point.  On  the  contrary,  using  the  tuning  strategy  of  Case 
II,  sharp  peaks  are  observed  as  far  as  the  dynamic  response  of  oxy¬ 
gen  excess  ratio  is  concerned,  but  at  the  same  time  the  dynamic 
behaviour  of  net  power  is  improved  (almost  instant  responses  to 
set-point  changes).  However,  there  is  one  exception  which  is  the 
net  power  response  when  the  set-point  is  set  to  its  highest  value. 


For  this  particular  set-point  step  change,  the  net  power  response  is 
as  slow  as  in  Case  I.  This  is  due  to  the  modeling  error  which  becomes 
significant  when  the  system  operates  in  a  region  which  is  not  close 
to  the  operating  range  that  was  used  for  obtaining  the  FIR  models. 

5.3.  Design  of  a  feedforward  control  scheme 

For  comparison  purposes,  a  feedforward  controller  was  also 
designed  based  on  the  results  of  Section  4.  In  particular,  the  infor¬ 
mation  contained  in  Fig.  11  was  used  to  develop  look-up  tables, 
from  which  the  optimal  values  of  the  manipulated  variables  are 
calculated,  as  functions  of  the  net  power  demand.  The  feedforward 
controller  was  applied  on  the  same  sequence  of  net  power  set-point 
changes  that  was  used  to  test  the  DMC  methodology  (Table  5).  The 
responses  of  the  controlled  variables  are  depicted  in  Figs.  17  and  18. 
It  is  clear  that  the  feedforward  controller  is  not  able  to  achieve  zero 
offset.  This  is  due  to  the  fact  that,  contrary  to  the  DMC  methodol¬ 
ogy,  the  feedforward  controller  lacks  integral  action.  Additionally, 
the  dynamic  behaviour  of  Aq2  contains  sharp  and  intense  changes. 
These  may  cause  the  abrupt  shut-down  of  the  system  in  a  real  appli¬ 
cation.  They  also  pose  dynamic  stress  onto  the  FC  membranes  due  to 
pressure  oscillations  and  possible  oxygen  starvation  [19,24],  which 
reduces  lifetime  of  the  membranes. 

5.4.  Testing  the  robustness  of  the  control  strategies 

In  this  subsection,  robustness  of  the  proposed  MPC  methodol¬ 
ogy  is  tested  on  a  typical  defection  case,  that  happens  during  the 
operation  of  fuel-cell  systems.  In  particular,  it  is  assumed  that  21 
of  the  381  cell  units  that  constitute  the  FC  stack  are  out  of  opera¬ 
tion,  due  to  membrane  cracking.  The  DMC  methodology  tuned  with 
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Fig.  15.  Response  of  net  power  (DMC,  Case  II). 


Fig.  17.  Response  of  net  power  (feedforward  controller). 
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Fig.  18.  Response  of  oxygen  excess  ratio  (feedforward  controller). 
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Fig.  21.  Response  of  oxygen  excess  ratio  in  robustness  test  (feedforward  controller). 
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Fig.  19.  Response  of  net  power  in  robustness  test  (feedforward  controller). 


Fig.  22.  Response  of  oxygen  excess  ratio  in  robustness  test  (DMC,  Case  I). 


the  Case  I  strategy  and  the  feedforward  controller  were  tested  on 
the  disturbed  system  and  the  responses  concerning  the  controlled 
variables  are  shown  in  Figs.  19-22.  Robustness  of  the  proposed 
method  is  illustrated  by  noticing  that  no  significant  modifications 
are  observed  in  the  dynamic  responses  of  the  controlled  vari¬ 
ables,  compared  to  the  case  where  all  cell  unit  operate  normally 
(Figs.  13  and  14),  and,  the  method  still  produces  zero  steady-state 
effort.  The  dynamic  behaviour  of  the  feedforward  controller  is  simi¬ 
lar  to  the  one  presented  in  Figs.  17  and  18,  but  the  steady  state  offset 
is  further  increased. 


Time  (s) 


Fig.  20.  Response  of  net  power  in  robustness  test  (DMC,  Case  I). 


6.  Conclusions 

In  this  work  a  computational  tool  for  the  optimal  control  of  FC 
systems  was  developed  by  integrating  and  combining  a  number 
of  modeling,  optimization  and  control  technologies.  The  produced 
methodology  is  built  upon  four  steps  that  may  be  summarized  as 
follows: 

•  First,  via  a  simulation  model,  a  database  is  generated  which 
includes  values  of  the  input  variables  (stack  current,  Jst,  compres¬ 
sor  voltage,  Ucm )  and  the  corresponding  values  of  output  variables 
(net  output  power,  Pnet ,  ratio  of  incoming  oxygen  mass  flow  to  the 
oxygen  that  reacts,  \<j2 ); 

•  Second  step  is  the  development  of  mathematical  relationships 
between  the  input  and  output  variables.  In  particular,  the  RBF 
neural  network  architecture  was  utilized  for  producing  the  input- 
output  “meta-models”; 

•  In  the  third  step,  an  NLP  optimization  problem  was  formulated 
that  takes  into  account  the  constraints  and  limitations  of  the  sys¬ 
tem.  The  objective  function  to  be  minimized  is  stack  current  for  a 
given  value  of  power  demand.  Based  on  the  formulation  and  solu¬ 
tion  of  the  NLP  problem,  a  look-up  table  is  developed,  containing 
the  optimal  values  of  the  system  variables,  for  any  possible  value 
of  power  demand  in  the  operating  range  of  the  system; 

•  Finally,  a  methodology  based  on  the  principles  of  MPC  is  designed 
for  the  optimal  control  of  the  system  response  to  successive 
sep-point  changes  of  power  demand.  The  optimal  Xq2 -value 
(determined  by  the  solution  of  the  optimization  problem  as 
described  above)  is  used  as  an  additional  control  variable.  The 
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efficiency  of  the  produced  MPC  system  is  illustrated  in  a  number 
of  simulations,  which  show  that  a  successful  dynamic  closed-loop 
behaviour  can  be  achieved,  while  at  the  same  time  the  consump¬ 
tion  of  hydrogen  is  minimized.  It  is  illustrated  that  by  appropriate 
tuning,  the  proposed  control  strategy  can  produce  zero  steady- 
state  offset,  without  the  risk  of  oxygen  depletion  (oxygen  excess 
ratio  does  not  fall  below  the  value  of  1.8  which  is,  in  fact,  much 
higher  than  the  safety  limit  of  1.5  [3]). 
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